Self-energy-functional approach: Analytical results and the Mott-Hubbard transition 
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The self-energy-functional approach proposed recently is applied to the single-band Hubbard 
model at half-filling to study the Mott-Hubbard metal-insulator transition within the most simple 
but non-trivial approximation. This leads to a mean-field approach which is interesting conceptually: 
Trial self-energies from a two-site single- impurity Anderson model are used to evaluate an exact 
and general variational principle. While this restriction of the domain of the functional represents 
a strong approximation, the approach is still thermodynamically consistent by construction and 
represents a conceptual improvement of the "linearized DMFT" which has been suggested previously 
as a handy approach to study the critical regime close to the transition. It turns out that the two-site 
approximation is able to reproduce the complete (zero and finite-temperature) phase diagram for the 
Mott transition. For the critical point at T = 0, the entire calculation can be done analytically. This 
calculation elucidates different general aspects of the self-energy-functional theory. Furthermore, it 
is shown how to deal with a number of technical difficulties which appear when the self-energy 
functional is evaluated in practice. 

PACS numbers: 71.10.-w, 71.15.-m, 71.30.+h 



I. INTRODUCTION 

The correlation-driven transition from a paramag- 
netic metal to a paramagnetic insulator (Mott-Hubbard 
transitioniSiii) is one of the most interesting problems 
in condensed-matter physics. As a prime example for a 
quantum-phase transition, the Mott-Hubbard transition 
is important from the physical point of view but also 
for the development and test of general theoretical meth- 
ods to treat correlated electron systems. The minimum 
model required to study the Mott-Hubbard transition 
is the single-band Hubbard modelA^ Inherent to this 
model is the competition between the electrons' kinetic 
energy which tends to dclocalize the electrons and fa- 
vors a metallic state and the on-site Coulomb interaction 
which tends to localize the electrons to avoid double oc- 
cupancies and thereby favors an insulating state at half 
filling. Except for the one-dimensional case^ however, 
exact results with regard to the nature of the transition 
and the critical interaction strength U c are not available 
- even for this highly simplified model system. A direct 
numerical solution using exact-diagonalization or quan- 
tum Monte-Carlo methods^ suffers from the difficulty to 
access the thermodynamic limit or the low-temperature, 
low-energy regimes. 

Considerable progress has been made in recent years 
due to the development of the dynamical mean-field 
theory (DMFT)&a2iii which focuses on the opposite 
limit of infinite spatial dimensions D = oo piiiiiii^i 
Within the DMFT the problem is simplified by map- 
ping the original lattice model onto an impurity model 
the parameters of which must be determined by a self- 
consistency condition. Different techniques to solve the 
effective impurity model have been employed to study 
the Mott transition within the DMFT, iterative per- 
turbation theory ; 9 i 16 i 17 exact diagonalization ) 18 i 19 i 20 i 21 
renormalization-group methods an H quantum 

Monte-Carloi 25 i 26 i 27 i 28 i 29 One of the most important 



characteristic of the transition is the value of the critical 
interaction strength U c at zero temperature. Roughly, 
the different techniques to solve the mean-field equations 
predict U c /W f» 1 — 1.5 where W is the width of the free 
density of states. 

Recently, a self-cnergy-functional approach (SEA) has 
been put forward^ The SEA is a general variational 
approach to correlated lattice models where the grand 
potential fl is considered as a functional of the self- 
energy S. As this functional is constructed from an 
infinite series of renormalized skeleton diagrams, it is 
not known in an explicit form and the variational prin- 
ciple = cannot be exploited directly. Usually, 
one replaces the exact but unknown functional with an 
explicitly known but approximate one - this is essen- 
tially the standard diagrammatic approach^ which leads 
to weak-coupling approximations in the end^i Opposed 
to this weak-coupling perturbational approach, the func- 
tional dependence f2[S] is not approximated at all in the 
SFA. The key observation is that the functional, though 
unknown explicitly, can be evaluated on a restricted do- 
main of trial self-energies S. The evaluation of the func- 
tional is exact, while the approximation is due to the 
fact that the self-energy in the variational principle is 
no longer considered as arbitrary. In this way, depend- 
ing on the choice for the space S, some well-known but 
also some novel approximations can be realized. Here 
we are interested in the single-band Hubbard model with 
Hamiltonian H. As argued in Ref. |3(] a useful trial self- 
energy has to be constructed as the exact self-energy of 
a different model ( "reference system" ) with Hamiltonian 
H' . The reference system can be chosen arbitrarily - it 
must, however, share the same interaction part with the 
original model H. The variational parameters at one's 
disposal are therefore the one-particle parameters of the 
reference system t' . The trial self-energy is parameter- 
ized as S = X(i'), and the variational principle reads 
dfl[E(t')]/dt' = 0. To provide trial self-energies is the 
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only purpose of the reference system H' . Whenever one 
is able to compute X for the reference system H' , an 
exact evaluation of 0[S(t')] is possible. 

Choosing H' to be a system of decoupled sites, yields 
a Hubbard-I-type approximation. An improved approx- 
imation is obtained when H' consists of decoupled clus- 
ters with a finite number of sites N c > 1 per cluster as 
has been considered in Refs. I3l3l This approach not 
only recovers the so-called cluster-perturbation theory 
(CPT)^ 4 * 3 ^ 3 ^ but also gives a variational improvement 
(V-CPT) which e.g. allows to describe phases with spon- 
taneously broken symmetry^ Another possibility is to 
take N c = 1, which implies the trial self-energy to be 
local, but to include a coupling to a number of rib addi- 
tional uncorrelated ("bath") sites. In this case the refer- 
ence system consists of a decoupled set of single-impurity 
Anderson models (SIAM) with n s = 1 + rib sites each. 
As has been shown in Ref. l3fl this approach not only 
recovers the DMFT (namely in the limit rib —* 00) but 
also provides a new variant of the exact-diagonalization 
approach, namely for any finite rib- As compared to pre- 
vious DMFT-exact-diagonalization approaches ^ItLZiLZL 
the construction gives a thermodynamically consistent 
approximation even for small rib- More complicated ref- 
erence systems may be taken for the construction of con- 
sistent approximations, for example a system of decou- 
pled clusters of size N c where each site in the cluster is 
coupled to rib additional bath sites. It has been shownS 
that in the limit n b ->• 00 the cellular DMFT (C-DMFT) 
is obtained^ while approximations with finite rib repre- 
sent cluster approximations which "interpolate" between 
the CPT (n b = 0) and the C-DMFT (n b = 00). 

In the present paper the Mott transition is studied 
within the most simple but non-trivial approximation: 
The variational principle <50[S] = is exploited using 
a local trial self-energy from a reference system with 
N c = 1 and a single additional bath site only, rib = 1. As 
n s = 1 + rib = 2 this approximation will be referred to as 
the two-site dynamical impurity approximation (ri s = 2- 
DIA) in the following. 

The paper is organized as follows: A brief review of the 
sclf-cncrgy-functional approach will be given in the next 
section m The general aspects of the evaluation of the 
SFA are discussed in Sec. IIIII while Sec. IIVI focuses on lo- 
cal approximations (N c = 1, n s arbitrary) in particular. 
In Sec. the further specialization to the case n s = 2 
is considered. The n s = 2 dynamical-impurity approx- 
imation is motivated (i) by the fact that at the critical 
point for the Mott transition the entire calculation can 
be done analytically, (ii) by the conceptual simplicity of 
the approach which rests on a single approximation only 
and (iii) by making contact with a linearized DMFT (L- 
DMFT) 38 i 39 i 40 i 41 i 42 i 43 developed previously. This is dis- 
cussed in detail in Sec. M while Sec. lVII then presents the 
analytical calculation for the critical regime. The results 
are discussed in Sec. IVIII The complete phase diagram 
for T = and finite temperatures is addressed in Sec. 
IVIIII and the conclusions are given in Sec. IIXI 



II. SELF-ENERGY-FUNCTIONAL THEORY 

Useful information on correlated electron systems can 
be gained by exact-diagonalization or quantum Monte- 
Carlo methods applied to a lattice of finite size£ This 
approach, however, suffers from the difficulty to access 
the thermodynamic limit and is therefore of limited use 
to describe phases with long-order and phase transitions. 
On the other hand, using an embedding approximation, 
one can directly work in the thermodynamic limit and de- 
scribe phase transitions while the actual numerical treat- 
ment has to be done a system of finite size only. In the 
context of an embedding technique we have to distinguish 
between the original model of infinite size H and an (e.g. 
spatially) truncated reference system H' . The reference 
system must not necessarily be finite but it may consist 
of an infinite number of decoupled subsystems with a fi- 
nite number of degrees of freedom each. In any case, H' 
must be exactly solvable. 

If one is interested not only in the equilibrium thermo- 
dynamics but also in the elementary one-particle excita- 
tions, an embedding technique should focus on a dynami- 
cal quantity, such as the frequency-dependent self-energy 
S. The knowledge of S then allows to derive a thermody- 
namic potential O as well as different static and dynamic 
quantities via general relations. The main steps are the 
following: (i) Truncate the original model H to obtain 
a simpler model H' which is tractable numerically, (ii) 
Calculate the self-energy X' of the reference system, (iii) 
Use S = £' as an approximation for the self-energy of 
H and determine the grand potential f2 as well as fur- 
ther quantities of interest, (iv) To get a self-consistent 
scheme, optimize the parameters of W by a feedback 
from the approximate solution at hand. Ideally, the last 
step should be based on a general variational principle. 
This is exactly the strategy of the self-encrgy-functional 
approach (SFA) proposed recently^ A brief review of 
the essentials of this approach is given in the following. 

Consider a system of fermions on an infinite lattice 
with on-site Coulomb interaction at temperature T and 
chemical potential fi. Its Hamiltonian H = Ho(t) + 
Hi(U) consists of a one-particle part which depends on a 
set of hopping parameters t and an interaction part with 
Coulomb-interaction parameters U: 

H = ^ tapclcp + - ^ UapSyC^pCyCS ■ (1) 
a/3 a0~f5 

The grand potential f2 can be obtained from the station- 
ary point of a self-energy functional 

fi t [£] = Trln(-(Go 1 - S)" 1 ) + F[£] (2) 

as has been discussed in Ref. |23 Here the subscript t 
indicates the parametric dependence of the functional 
on the hopping. This dependence is exclusively due to 
Go = + // — *), the free Green's function of H . Fur- 
ther, F[S] is the Legendre transform of the Luttinger- 
Ward functional 3>[G]. As the latter is constructed as 
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(a) (b) 

FIG. 1: Schematic representation of the single-band Hubbard 
model H (a) and a possible reference system H' (b). H' is 
a set of decoupled single-impurity Anderson models with one 
correlated (U > 0) impurity site and a number of n s — 1 
uncorrelated (U = 0) bath sites each. In the figure n s = 4. 
Note that the interaction part is the same for (a) and (b). 
Variational parameters are the one-particle parameters of H' . 

an infinite series of renormalized skeleton diagrams^ the 
self-energy functional is not known explicitely. Neverthe- 
less, the exact evaluation of fit [S] and the determination 
of the stationary point is possible^ on a restricted space 
S of trial self-energies E(i') e S. Due to this restriction 
the procedure becomes an approximation. 

Generally, the space S consists of t 1 representable self- 
energies. X is termed t' representable if there are hopping 
parameters t' such that £ = S(t') is the exact self-energy 
of the model H' = H (t') + Hi(U) ("reference system"). 
Note that both the original system H and the reference 
system H' must share the same interaction part. For any 
X parameterized as S(t') we then havei^ 

fit [£(*')] = fi'(t') 

+ TVlnHGoCt)- 1 -^'))- 1 ) 
- TrlnC-CGoCt'r 1 -^'))- 1 ), (3) 

where fi'(t'), G' = G {t') = l/(u +n~ t'), and S(t') 
are the grand potential, the free Green's function and 
the self-energy of the reference system H ' while Go is 
the free Green's function of H . For a proper choice of 
t' , a (numerically) exact computation of these quantities 
is possible. Hence, the self-energy functional J2J) can be 
evaluated exactly for this X = E(£'). A certain approxi- 
mation is characterized by a choice for S. As S is param- 
eterized by t', this means to specify a space of variational 
parameters t' . Any choice will lead by construction to a 
non-perturbativc approach which is thcrmodynamically 
consistent as an explicit expression for a thermodynam- 
ical potential is provided. It turns out that a stationary 
point of the self-energy functional is a saddle point in gen- 
eral. As in different standard variational methods, such 
as in the time-dependent density-functional approach^ 
in the Green's-function approach^! and also in a recently 
considered variant^ this implies that there is no strict 
upper bound for the grand potential. For a further dis- 
cussion of the general concepts of the SFA see Ref. |3?1 

So far the discussion is completely general. In Sec. 
IIVI we will consider H to be the Hubbard model and 
H' to be a system of decoupled single-impurity Ander- 
son models (SIAM). Each SIAM consists of n s sites, one 



correlated site (with U > 0) and n s — 1 uncorrelated 
"bath" sites (U = 0). This is illustrated by Fig. [T] for 
n s = 4. Note that for any choice of n s , the original 
system and the reference system share the same inter- 
action part - as required by the general theory. The 
one-particle parameters of H' are the variational param- 
eters, i.e. the one-particle energies of the original sites and 
the bath sites and the hop ping ( "hybridization" ) between 
them. As noted in Ref. |3fJ, the dynamical mean-field 
theory (DMFT) is recovered in the limit n s — + oo. For 
n s < oo one obtains a new variant of the DMFT-exact- 
diagonalization approach>i2ii2*22i2i For a finite number 
of bath sites, the DMFT self-consistency condition can- 
not be strictly satisfied. In the DMFT-ED method one 
therefore has to introduce a certain measure which al- 
lows to minimize the error due to the discretization of 
the bath. The conceptual advantage of the SFA con- 
sists in the fact that this measure is replaced by a varia- 
tional procedure which is based on a physical variational 
principle. As shown in Ref. l30t a very good quantita- 
tive agreement with results from full DMFT calculations 
can be achieved for the quasi-particle weight even with 
n s = 4. With n s ~ 2 & much simpler approach is con- 
sidered here (Sees. |V] and IVIf) which, however, is still 
thermodynamically consistent and allows for simple and 
systematic investigations of the Mott transition. 



III. EVALUATION OF THE SELF-ENERGY 
FUNCTIONAL 

The evaluation of the self-energy-functional theory can 
be done by solving the Euler equation^ corresponding 
to the variational principle. While such an approach is 
possible in principle, it appears complicated as dynam- 
ical two-particle quantities of the reference system are 
required. An attractive alternative consists in the direct 
calculation of the grand potential along Eq. ©. The 
numerical computation of fit[S(t')] for a given set of 
one-particle parameters t' is straightforward for a ref- 
erence system H' of finite size. There are, however, a 
few technical difficulties which appear in the practical 
calculation and which shall be discussed in the following. 
The problem of finding a stationary point of the function 
t' — * fit[S(t')] is not addressed here as this is a standard 
numerical problem very similar to the problem of find- 
ing a minimum of a real single- valued function of several 
arguments. 

fit[S(t')] consists of three parts as given by Eq. 
The grand potential of the reference system can be cal- 
culated as fi'(i') = -Tlntr'exp(-(if' - ^N')/T) = 
- Tln E,„ cx P(-(^m - t+N'J/T) from the many-body 
cigenencrgies E' m — ^N' m of H' — fiN' where N' is the 
total particle number operator. Direct numerical diago- 
nalization or (at T = 0) the Lanczos technique^ may be 
used. 

Next, the second term on the r.h.s. of Eq. (0 is dis- 
cussed. In the following the dependence of X on t' will 
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be suppressed for convenience and its dependence on lo 
is made explicit in the notations. The diagonalization 
of H' — fxN' yields (via the Lehmann representation) 
the Green's function G'(lo) and the free Green's function 
G' (ui) of the reference system. The self-energy is then 
obtained from the Dyson equation of the reference sys- 
tem = G' Q 1 (lo) — G' 1 (lo). Using the self-energy 
one obtains the (approximate) Green's function of 
the original model via G(lo) = (Gq 1 ^) — S(oj)) -1 . The 
remaining task is to calculate Tr ln(— G) for a lattice con- 
sisting of a finite number of sites L. The thermodynamic 
limit L — > oo is performed in the end. Translational 
symmetry is not necessarily required. 



It is important to note that G(u>) is causal, i.e. G{u> + 
i0 + ) = Gr(w) — iG\{uj) with Gr(w), G\{lo) Hermitian 
and Gj{u>) positive definite for any real lo (0 + is a posi- 
tive infinitesimal). The causality of the Green's function 
G(iS) is ensured by the causality of S(w) and Gq(u>), see 
Appendix The latter are causal as these are exact 
quantities. 



Let u) m be the (real, first-order) poles of G{io). For 
lo — > uj m we have G(lo) — > R m /(uj — Lo m ) where due 
to the causality of G{lo) the matrix R m is positive defi- 
nite. Therefore, the frequency-dependent diagonalization 
of G(uj) = U(Lo)g(Lo)U(Loy with unitary U(lo) (for real 
to) yields a diagonal Green's function g(u>) with elements 
gk(j-o) that are real for real lo and have first-order poles 
at to = to m with positive residues. 

As the Green's function can be written as G(to) = 
l/(to + — t — S(cj)), the unitary transformation U(uj) 
also diagonalizes the real and symmetric matrix t+E(a;), 
i.e. t + = U(to)rf(to)U(to)^ and g k {to) = l/(w + fi — 
77fc(w)). The rjk{io) are real for real to and have first- 
order poles at to = £ n with being the poles of the 
self-energy For to — > £„ we have S(cj) — > S n /(uj — £„) 
with positive definite S^. Consequently, the residues of 
i]k(to) at w = Cn are positive. 

One can write g- fc (w) = Em^,™/( w_w "')+ft( w ) and 
?7fe(w) = ^ n Sk,n/(u ~ Cn) + ?7fc(w) with Rk, m , Sk,n > 
and where g k {to) and rjk(to) are analytical in the en- 
tire to plane. As G(lo) ~ l/u> and ~ const, 
for cj — > oo, one has g k (Lo) = and rj k {to) = const, 
and real. Consequently, — (l/ir)Imgk(to + i0 + ) > and 
-(l/7r)Imr?fe(w + £0+) > 0. This will be used in the fol- 
lowing. 

The trace "Tr" in Eq. © consists of a sum T ^ w over 
the fermionic Matsubara frequencies ito = i(2n + 1)ttT 
(n integer) and a trace "tr" with respect to the quantum 
numbers a; see Eq. Q . The convergence of the frequency 
sum is ensured by the usual factor cxp(ia>0 + ) from the 



diagram rules. The calculation then proceeds as follows: 

-1 



TVe iu0 tr In - 



ito + fj, — t — S(icj) 



(a) -1 



k c 

= —Y f°° du f{to) Im ln(-g k (Lu + »0+)) 

( c ) Z" 00 

7. J — OO 



(d) 
(J 

with 



y r^T in(i+ e -^) + ; -^M) 

-2L ^ T ln(l + exp(-w m /T)) - ifc (4) 



JJ s = -^Tln(l + exp(-C„/T)). (5) 



In Eq. J2J In denotes the principal branch of the loga- 
rithm, f(to) = l/(exp(w/T) + 1) is the Fermi function, 
and C is a contour in the complex to plane enclosing 
the first-order poles of the Fermi function in counter- 
clockwise direction. In step (a) the transformation U(to) 
is performed under the trace. Convergence of the inte- 
gral for to — > ±oo is ensured by the Fermi function and 
by the factor e w0+ , respectively Step (b) results from 
analytical continuation to real frequencies. In step (c) 
— (l/7r)Im(7fc(u; + i0 + ) > has been used (see Appendix 
IB1 . At this point the causality of the Green's function is 
essential, as discussed above. In step (d) the Fermi func- 
tion is written as a derivative with respect to lo and in- 
tegration by parts is performed. Step (e) uses the results 
of Appendix [U] for the derivative of the step function. As 
the different diagonal elements of the Green's function, 
g k (tj), have the same set of poles and zeros^ the sum 
over k becomes trivial and yields a factor 2L only where 
L is the dimension of the hopping matrix, i.e. the num- 
ber of orbitals. The factor 2 accounts for the two spin 
directions. The contribution from the poles of the self- 
energy (Appendix Q is denoted by Apart from this 
correction term, Tr ln(— (Go(i) -1 — S(t')) -1 ) turns out 
to be the grand potential of a system of non-interacting 
quasi-particles with unit weight and energies given by the 
poles of G(lo) = (Gq 1 (to) - S(w)) _1 . 

Consider now the third term on the r.h.s. of Eq. 
A calculation completely analogous to Eq. results in: 

T^e iw0+ tr' ln(-G'M) 

UJ 

= -2i^Tln(l + exp(- W ; i /T))-i? s . (6) 

m 

Here Lo' m are the poles of G'. Again, the first term in 
Eq. © is the grand potential of a non-interacting sys- 
tem of fcrmions with onc-particlc energies given by the 



5 



poles of the Green's function G . The same holds for the 
second term, but with energies given by the poles of the 
sclf-encrgy. By construction, the self-energy is the same 
for both, the original system and the reference system. 
Hence, the same correction term i?s appears in Eqs. 
and © and cancels out in Eq. 0. 

Note that a pole of G(w) at lu = Lu rn with residue 
Rrn — * implies a pole of G'(u>) at the same frequency 
w,„ (with residue R' m — > 0). Hence, contributions due 
to poles with vanishing residues cancel out in Eq. Jjyl. 
The reason is the following: Suppose that gk(^) = 
Rk,m/ — w m ) for lu close to tu m with residue Rk, m —> 0. 
For the diagonal elements of t + this implies that 

r) k (uj) = (l/i4 >m )(w-w m ) near w m . A zero of rj k (ui) with 
infinite positive coefficient l/Rk,m must be due to S(w). 
Therefore, for the diagonal elements of t' + this 
means that r)' k {io) = (1/ ' Rk, m ){to — w m ). Consequently 
g' k {uj) = l/(u) + p, - rj' k {uj)) = Rk, m / (w - oj m ) for uj close 
to uj m with R k _ m — > 0. The argument also shows that 
although the residues do not appear in Eqs. |J3J and (JSJ 
explicitly, one can state that poles with small residues 
will give a small contribution in Eq. J3J| . 

Based on the causality of the respective Green's func- 
tions, an efficient algorithm can be set up to find ui m 
and ui' m numerically which are then needed in Eqs. Q 
and JfjJ). The poles of G'(ui) are directly obtained from 
the diagonalization of the reference system. The prob- 
lem consists in finding the poles of G(lu). Allowing for 
poles with vanishing (very small) residue, one can as- 
sume the function g' k {u)) for fixed but arbitrary k to dis- 
play all the poles of G'(uj). Since g k (uj) is of the form 

0k M = Em R'k,mJ ( W ~ ^m) With R' k>m > 0, it is HlOnO- 

tonically decreasing. Hence, there is exactly one zero of 
g' k {w) located in the interval between two adjacent poles 
u>' m and (jj' m+ i- As g' k (u>) is monotonous, the zero £ m can 
easily be found numerically by an iterative bisection pro- 
cedure. The zeros of g k (io) are the poles of r] k (ui) and 
the poles of r]' k (ui) are the same as the poles of r) k {io). 
Now, since g k (u) = l/(w + fJ, - %(w)), #fc(o>) and g' k (co) 
must have the same set of zeros. The function g k {uj) is 
monotonically decreasing. Therefore, the poles of g k {u>) 
can be found between the ( m by using the same iterative 
bisection procedure once again. 



This is the DMFT which is recovered for an infinite num- 
ber of uncorrelated bath sites (per original correlated 
site), i.e. for n s — 1 — * oo. On the other hand, n s = 1 
corresponds to a Hubbard-l-type approximation. Here 
it will be shown that, for arbitrary n s , the evaluation of 
the self-energy functional reduces to a one-dimensional 
integration only. This is an important simplification for 
any practical numerical (or even analytical) calculations. 

In case that the self-energy is local it is advantageous 
to start from step (c) in Eq. Q). Assuming translational 
symmetry, the matrix t + H(uj) is diagonalized by Fourier 
transformation to reciprocal space. Its eigenvalues rj k (uj) 
are given by rj k (uj) = e(k) + E(w) where k = (k,a) 
and e(fc) is the tight-binding Bloch dispersion. The self- 
energy is taken to be spin-independent and independent 
of the site index, i.e. a paramagnetic homogeneous phase 
is assumed for simplicity. The second term on the r.h.s. 
of Eq. 10 then becomes: 



TrlnHGoCt)- 1 -^'))" 1 ' 

(c) 



= -Yj J duj /H 9 ( u + m - £ ( fc ) - S H) 

ka 

2L / dw f(ui) / dz po(z) 8 (uj + /i — z — S(w)) 



2L I dw f(uj) I dz R (z) —6 (cj + /a - z - E(w)) 



2L / duj f{u>) R {uj + £(w)) 



(8) 



Here L is the number of lattice sites, the factor 2 stems 
for the spin summation, po(z) = L~ 1 ^2 k 5(z — e(fc)) 
is the non-interacting density of states, and Ro{z) — 
J^^dz'poiz') is an antiderivative. Note that there is 
no "correction term" as the derivative in Eq. (JSJ is 
with respect to z. 

The reference system H' is a set of decoupled single- 
impurity Anderson models with one correlated ("impu- 
rity") site and n s — 1 bath sites each. The Hamiltonian 
is H' = with 



IV. LOCAL APPROXIMATIONS 

In the following we consider H to be the single-band 
Hubbard model: 



H 



ijcr 



2 ^ 



(7) 



u 



H'(i) = E e ^ + yE 



E £k4 k<T a ika + (VkctrOik* + H.c 



;k=2 



(9) 



For the reference system H' shown in Fig. ^ t ne self" 
energy is local: Sy(w) = 5yE(o>). Clearly, the approx- 
imation is the better the more degrees of freedom are 
included in H'. The optimal local approximation is ob- 
tained with the most flexible (but local) trial sclf-encrgy. 



The parameters e±, e k , and V k for fc = 2, ...,n s are the 
on-site energies of the impurity site and of the bath sites, 
and the hybridization between them, respectively These 
are the variational parameters of the theory 
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The third term on the r.h.s. of Eq. (J3J) reads: 
TrlnHG^'r 1 -^'))- 1 ) 



fc=2 



(L/m e(GiM) + ^e(Gl.H) 

00 V fc=2 



(10) 



where the first equation is derived in Appendix IdI and 

1 



G'M = 



uj + fi — t\ — A(ll>) — £(w) 
is the impurity Green's function, 

1 



u + /i - e k 

is the (free) conduction-band Green's function, and 

A(-) = E 



k=2 



cj + fx- e k 



(11) 



(12) 



(13) 



is the hybridization function. The final expressions © 
and I|1U|) involve one-dimensional integrations only and 
can therefore be calculated numerically without serious 
problems. 



V. TWO-SITE DYNAMICAL IMPURITY 
APPROXIMATION 

In the following we will focus on the case n s = 2, i.e. on 
the two-site dynamical-impurity approximation (n s = 2- 
DIA). There are different intentions which are followed 
up: 

(i) For any approximation within the context of the 
SFA, one has to compute the self-energy of the refer- 
ence system H' . As the interaction part is the same for 
both, H and H' , this still constitutes a non-trivial many- 
body problem which can only be treated by numerical 
means in most cases. The reference system characterized 
by N c = 1 and n s = 2, represents an exception: For a 
special point in the space of model parameters (zero tem- 
perature, half- filling and U = U c , the critical interaction 
for the Mott transition), the entire calculation can be 
done analytically. This not only includes the diagonal- 
ization of the reference system which actually is a simple 
dimer model but also and more important here, the ex- 
act evaluation of the self-energy functional for the trial 
self-energies considered and the subsequent variational 
optimization. Therefore, the study of the n s = 2-DIA is 
ideally suited to elucidate different technical points which 
are relevant for any A^ c and n s and which must be con- 
sidered carefully. 

(ii) The n s = 2-DIA must be considered as inferior 
when compared to approximations with higher n s and 



when compared to n s = 00 (the DMFT) in particu- 
lar. On the other hand, one has to keep in mind that 
the DMFT must always be supplemented by an addi- 
tional (numerical) method to solve the mean-field equa- 
tions which necessarily involves additional approxima- 
tions. Even if the additional approximations can be con- 
trolled, DMFT results always depend on the accuracy of 
the numerical method employed. As concerns the n s = 2- 
DIA, there is no such difficulty: The theory rests on a 
single approximation only, namely the restriction of the 
space S to self-energies representable by the two-site ref- 
erence system - the rest of the calculation is rigorous. It 
is this conceptual simplicity which makes the approxima- 
tion attractive. 

(iii) That an approach referring to an n s = 2-sitc SIAM 
is able to give reasonable results has been shown be- 
forehand by the linearized DMFT (L-DMFT) , 38 i 39 i 40 i 41 i 42 
The L-DMFT is a well-motivated but ad-hoc simplifica- 
tion of the full DMFT and maps the Hubbard model self- 
consistcntly onto the n s = 2-site SIAM just at the criti- 
cal point for the Mott transition. The L-DMFT can be 
also be considered to represent the lowest-order realiza- 
tion of a more general projective self-consistent method 
(PSCM)^ 2 . As compared to the full DMFT, the linearized 
theory yields surprisingly good estimates for the critical 
interaction U c in the single-band model, on translation 
invariant lattices^ as well as for lattice geometries with 
reduced translational symmetry^ 9 . The approach can also 
be extended beyond the critical regime^ 3 - The main dis- 
advantage of the L-DMFT is that it not consistently de- 
rived from a thermodynamical potential. Another inten- 
tion of the present paper is therefore to suggest a two-site 
method that conceptually improves upon the L-DMFT 
in this respect. In fact, as the parameters of the effec- 
tive n s = 2-site impurity model are determined via a 
physically meaningful variational principle, the two-site 
approximation within the SFA should be regarded as an 
optimal two-site approach. The interesting question is, of 
course, whether or not this improvement of the method 
also implies improved results. To this end the results 
from the analytical evaluation of the n s = 2-DIA have to 
be compared with those of the L-DMFT and with avail- 
able numerical results for n s ~ 00 (full DMFT). This 
represents a good check of the practicability of the new 
method. As the n s = 2-site DIA still represents a very 
handy method, it can also be employed to investigate 
overall trends. It will be interesting to study the critical 
interaction U c for a variety of different geometries (i.e. 
for different free densities of states). 

Let us first consider a numerical evaluation of the the- 
ory. Fig. [5] shows the results of a numerical calcula- 
tion along the lines discussed above for the paramagnetic 
phase of the Hubbard model at half-filling and zero tem- 
perature. The free density of states (DOS) po(z) is taken 
to be semi-elliptical with a band width W = 4. The cal- 
culations are performed using the n s = 2-DIA, i.e. there 
is one bath site (per correlated site) only. Due to man- 
ifest particle-hole symmetry, two of the variational pa- 
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FIG. 2: The grand potential Q (per site) as calculated from 
Eq. © and Eqs. JHJ and 1101 for the reference system with 
n s = 2 as a function of the hybridization strength V (only 
the difference fl(V) — fi(V = 0) is plotted). Calculations 
for the Hubbard model at zero temperature T = and half- 
filling (/i = U/2). The free density of states is taken to be 
semi-elliptical with a band width W = 4. The interaction 
strength is varied from U — 5 to U = 6 as indicated. Arrows 
indicate stationary points of the self-energy functional. The 
Mott transition takes place at a critical interaction strength 
U c « 5.85. 



rameters, the on-site energies, are already fixed: ei = 
and 62 = fi = U/2. As a function of the remaining vari- 
ational parameter, the hybridization strength V = Vk=2, 
the grand potential fl(V) = f2[S(V)] shows two (non- 
equivalent) stationary points for U = 5 (see Fig. 0): 
a minimum at a finite V — 0.37 and a maximum at 
V = (as fi(y) = f2(— V), there is another minimum 
at V = —0.37 which can be ignored here). The two 
stationary points correspond to two physically different 
phases: For V > the interacting local density of states 
is finite at u> = while it vanishes for V = 0. So there is 
a metallic and an insulating phase coexisting. Due to the 
lower at the respective stationary point, the metallic 
phase is stable as compared to the insulating one. With 
increasing U the optimal V mo t and the energy difference 
l^(Knct) - fi(0)| decrease. For U = U c w 5.85 there 
is a metal-insulator (Mott-Hubbard) transition which is 
characterized by a coalescence of the stable metallic with 
the metastable insulating phase. For U > U c there is the 
insulating phase only. 

Qualitatively, this continuous transition is completely 
consistent with t he preform ed-gap scenaric«ii*2£ (how- 
ever, see also Refs. l47l48l49j) . For U < U c the self-energy 
is a two-pole function. This leads to a three-peak struc- 
ture in the interacting local density of states: Similar 
as in the full DMFT, there are two Hubbard "bands" 
separated by an energy of the order of U, and a quasi- 
particlc resonance at uj = 0. On approaching the critical 
interaction U c from below, the weight z of the resonance 
vanishes linearly z ~ (U c — U) leaving a finite gap for 



U > U c . As shown in Ref. l30l the quasi-particlc weight 
calculated from the self-energy at the respective optimal 
V = Vmet is in a very good quantitative agreement with 
results from full DMFT calculations in the whole range 
from U = to U = U c - for n s = 4 and even for n s = 2 
which is the case considered here. 



VI. "LINEARIZED" DYNAMICAL IMPURITY 
APPROXIMATION 

In the following we will concentrate on the critical 
regime U — * U c . It will be shown that the critical inter- 
action strength can be calculated analytically for n s = 2. 
The independent analytical result can be compared with 
the numerical one of the preceeding section. This repre- 
sents a strong test of the numerics. 

Consider the function Q(V) = fi[E(V)]. As fl{V) = 
fi(— V) there must be a stationary point of Q(V) at V = 
for any U. This implies that the linear term in an 
expansion around V = is missing, i.e.: 

n(V) = fi(0) + A ■ V 2 + C(y 4 ) . (14) 

The coefficient A depends on U. Assuming that the metal 
is stable against the insulator for U < U c , we must have 
A < for U < U c and A > for U > U c . This is 
a necessary condition for a continuous ("second-order") 
transition and consistent with the numerical results dis- 
played in Fig. |21 Therefore, the critical point for the Mott 
transition within the two-site model is characterized by 

A = 0. (15) 

The task is to calculate the three contributions to the 
grand potential following Eq. J3J), to expand in V up 
the the second-order term and to find the interaction 
strength satisfying the condition 11511 . 

Consider the grand potential of the reference system 
first. With ei = 0, e 2 = U/2 and fi = U/2 the ground 
state of the two-site system, Eq. @, lies is the invariant 
subspace with total particle number N' = 2. The ground- 
state energy E' is readily calculated: 

E' = hi- ^VipTmV 2 . (16) 

At T — the grand potential of the reference system is 
n' = LE' Q - Lfi(N'). Therefore, 



n'/L = - u -- -y/u* + Mv* = -^-^ r + OiV') . 

(17) 

Actually, this is the grand potential per site. 

We proceed with the third term on the r.h.s. of Eq. 
(jSJ). For the analytical calculation, it is convenient to 
start from Eq. ©: 

T^tr' ln(-G'M) = 2L^J r Q{-u' r ) - i? s (18) 

uj r 
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where it has been used that — T ln(l + cxp(— uj/T)) = 
w9(-w) for T = 0. The factor 2 is due to the spin de- 
generacy. i?s will cancel out later. The Green's function 
of the two-site model is easily calculated. There are four 
excitation energies, labeled by r, given by the four poles 
of the impurity Green's function at: 



±- (V^ 2 + 64V 2 ± ^U 2 + 16V 2 ) . (19) 



This yields: 



T^tr' \n(-G'(iuj))/L = -17- 



32V 2 



-0{V A )-Rv/L. 

(20) 

Finally, for the second term on the r.h.s. of Eq. we 
have: 

T^tr ln(-G(icj)) = 

/oo 
dz p (z) uj r (z)Q(-uj r (z)) - R s ■ (21) 
-00 

Here, Eq. has been used with k = (k,a), and the k 
sum has been replaced by an integration over z weighted 
by the free DOS po{z). For a given z = e(k), the quasi- 
particle energies are obtained as the poles of the lattice 
Green's function in reciprocal space, i.e. from the solu- 
tions u) = u) r (z) of the equation u) + p — z — = 0. 
To find the roots, the self-energy of the reference system 
is needed: 



EH 



u u 2 



1 



1 



(22) 



2 8 \oj-3V oj + 3V, 

This leads to a cubic equation: 

us 3 - zlu 2 - (9V 2 + U 2 /4)uj + 9zV 2 = . (23) 

The solutions for V = are easily obtained. For small 
V =/= we find one root near to = 0: 



(24) 



and another one near ui 
z 1 



-U/2: 



18z 18 z 2 + U 2 /2 

U 2 ~ + u 2 ^/F+W 



v 2 + o(v 4 ) . 



(25) 



Because of the step function in Eq. I|21|l . the third root 
near ui = U/2 is not needed here. This yields: 

T^tr \n{-G(iuj))/L = -R s /L + 

2 J°° dz p (z) (e(z)^ V 2 + Lo 2 (z)^j + 0(V 4 ) . 

(26) 



Inserting the results, Eqs. (|T7jl . l(2Tj|) . and into Eq. 
and using the symmetry po(z) = Pq(—z), we find: 

(24 72 f° 
JJ + Jp J dz po(z) z 



36 

U 2 



(27) 



Now, the condition l|15|) gives the critical U for the Mott 
transition: 

f° 3 f°° z 2 + U 2 /2 

U c = -3 dzp (z)z + - dzp {z)—= 2 c . 

J —OO J — OO \J Z ~\- U Q 

(28) 

This implicit analytical equation for U c is the final result. 
For an arbitrary free DOS no further simplifications are 
possible. 



VII. DISCUSSION 

It should be stressed once more that Eq. results 
from an exact variational principle simply by the restric- 
tion that the trial self-energies be rcprcsentable by the 
two-site reference system. The two-site model generally 
yields a two-pole self-energy which is the minimal re- 
quirement for a three-peak structure of the single-particle 
excitation spectrum. Therefore, it may also be stated 
that Eq. (|28|) gives the optimal result for a two-pole self- 
energy. 

Eq. (|28|) turns out to be more complicated as compared 
to the result of the linearized DMFTiiS 



[/(L-DMFT) = g 



dz po(z) z 2 



(29) 



Interestingly, the first term on the r.h.s. of Eq. I|28|l re- 
sembles Brinkman and Rice (Gutzwiller) result^ for the 
critical interaction: 



c/i BR ) 



-16 



dz po(z) z . 



(30) 



A rough approximation of the second term in 128fl us- 
ing (z 2 + U 2 /2)/y/z 2 + U 2 -> U c /2 then yields U c w 

3t/i BR V4- This reduction of the critical interaction as 
compared to the BR result is the dominating effect. 

For a semi-elliptical free DOS with variance A = 1 
(band width W = 4), one has Uc L DMFT - ) — 6 and 
£/ c (BR) = 6.7906. The numerical solution of Eq. (j2EJ is 
straightforward and yields U c = 5.8450. As expected 
this is fully consistent with the numerical determination 
of U c (Sec. [VJ. The result is close to the one predicted 
by the L-DMFT. Compared to numerical estimates from 
full DMFT calculations u^ BMFT) = 5.84 (NRG, Refill 
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FIG. 3: Critical interaction U c as a function of the inverse 
band width W -1 . Calculations for different (normalized) free 
densities of states po(z) as indicated (where —3 < z < 3 and 
< po(z) < f.5 in the insets). 



L-DMFT yields U c = 6 irrespective of the form of the 
DOS. 

The critical interaction from the two-site DIA is always 
close the L-DMFT result but considerably lower than the 
Gutzwiller value. The two-site DIA confirms the central 
prediction of the L-DMFT that it is the variance of the 
DOS that is crucial for the critical interaction. However, 
there is also a weak trend superimposed, namely a sys- 
tematic increase of U c with decreasing band width W 
(with the exception of the Gaussian DOS). This is the 
same trend that is also present in the Gutzwiller results. 
It would be interesting to see whether or not this trend 
is confirmed by full DMFT calculations. Comparing the 
full DMFT results for the Gaussian, for the rectangular 
and the semi-elliptic DOS shows the mentioned trend. 
However, the comparison of the two-site DIA and of the 
L-DMFT with the full DMFT for the available numerical 
data is not fully conclusive. 



VIII. FINITE TEMPERATURES 



TABLE I: Critical interaction U c as obtained from the SFA 
within the n s — 2 dynamical impurity approximation for 
different free densities of states. Each DOS is normalized, 
Po(z) = 1, symmetric, po(z) — p (—z), and has unit vari- 
ance, A = f™ x z 2 Po (z) = 1. A(z) = Q(W/2 + z)Q(W/2 - z) 
is a cutoff function. W is the band width. The critical inter- 
action from the Gutzwiller (Brinkman-Rice) approach— and 
the available full DMFT results are shown for comparison. 
The linearized DMFTS yields ry c ( L - DMFT ) = 6 f or any DOS 
with unit variance. 



Po(z) 



W U? 



BR 



E4 U 



Gaussian oo 5.5663 6.3831 

3A(z)(\z\ - \/T0) 2 /20Vl0 2V10 5.5284 6.3246 

triangle 2\/6 5.6696 6.5320 

semi-ellipse 4 5.8450 6.7906 

rectangle 2^3 5.9385 6.9282 

A(z)\z\/2 2V2 6.3554 7.5425 

9V3A{z)z 2 /10V5 2/573 6.4944 7.7460 



5.80° 



34°,5.i 
6.32 c 



S(z- l)/2 + <5(z + l)/2 



6.6686 8.0000 



C R. Bulla, private communication. 



and U, 



(DMFT) 



5.88 (PSCM, Ref.|22j), one can state that 



the prediction of the L-DMFT is improved. 

Results for different free densities of states are dis- 
played in Fig. |2| and Tab. [I] For a meaningful com- 
parison, each DOS has unit variance A = 1 where 
A 2 = J^ c oc dzp (z)z 2 . The band width W varies. The 
DOS with the smallest W (but A = 1) consists of two 
(5-peaks while W = 00 for a Gaussian DOS. Clearly, there 
is no true Mott transition in the former case (but also in 
the latter this is questionable). However, the inclusion 
of these extreme cases is instructive when studying the 
trend of U c as a function of W. Note that for A = 1 the 



So far only the zero-temperature limit has been consid- 
ered. The Mott transition at hnitc temperatures, how- 
ever, is particularly interesting as there is a comprehen- 
sive physical picture available with a comparatively com- 
plex phase diagram. This phase diagram in the U-T 
plane was first suggested by the iterative perturbation 
theoryAii The nature of the transition and the topology 
of the phase diagram have been established entirely us- 
ing analytic arguments^ and has been worked out quan- 
titatively using different numerical methodsi^iiSLSSiSi 
The critical regime is accessible by the projective self- 
consistent method^ 

As this phase diagram represents a valuable benchmark 
for any approximation, it is interesting to see whether 
or not it can be rederived within the most simple two- 
site DIA. The application of the theory for finite T is 
straightforward but can no longer be done analytically. 
As for the derivation of the T = results in the non- 
critical regime (see Sec. 0), calculations are performed 
along the lines of Sec. II VI 

It turns out that for finite temperatures the Mott tran- 
sition is predicted to be discontinuous. This is demon- 
strated in Fig. 21 which shows fl(V) for different T and 
fixed U = 5.2 (n s = 2). At low but finite T there are three 
stationary points (see arrows) corresponding to three dif- 
ferent phases of the system. The metallic phase has the 
largest value V = V me t- With increasing temperature 
Vmct decreases and £l(V me t) oc T 2 for low T. This gives 
a linear entropy S(T) = —d£l(T)/dT and a linear spe- 
cific heat C v = TdS(T)/dT = 7T cx z^T as it is char- 
acteristic for a Fermi liquid (p, = U/2 is fixed). The 
insulating phase has the smallest value V = Mns, and 
the grand potential Q,(V ins ) oc T for low T. For T -> 
the entropy approaches S — > L In 2 reflecting a 2 L -fold 
ground-state degeneracy of the insulator which is known 
to be an artifact of mean-field theory>ii The specific heat 
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FIG. 4: Grand potential per site at U = 5.2 for different T 
as a function of V . The arrows indicate the stationary points. 
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FIG. 5: Phase diagram for the Mott transition in the para- 
magnetic phase of the Hubbard model at half-filling. Calcu- 
lations for a semi-elliptical density of states with band width 
W = 4. Reference system H': two-site model. U C 2'- criti- 
cal interaction strength up to which there is a (metallic) so- 
lution smoothly connected with the metal at U = 0. U c i- 
critical interaction strength down to which there is an (insu- 
lating) solution smoothly connected with the Mott insulator 
for U — ► oo. U c : first-order transition line in the coexistence 
region, terminating at the second-order critical point at T c . 



C v vanishes exponentially for T — > 0. For U = 5.2 and 
low T, however, the insulating phase is metastable as 
compared to the metallic phase since 0(Vj ns ) > fi(^met)- 
Due to the different behavior at low T, there must be 
a temperature T C (U) where f2(V; ns ) = Q(V me t)- In fact, 
for T = 0.012 > T C (U = 5.2) the insulator is stable as 
compared to the metal. At T C (U) or, conversely, at a 
critical interaction U C (T), the n s = 2-DIA thus predicts 
a first-order transition with a discontinuous jump in the 
entropy. As can be seen in Fig. the metallic phase 
ceases to exist for still higher temperatures, and only the 
insulating phase is left. This is due to a coalescence of 
the metallic with a third phase at another critical tem- 
perature T c2 (U) (or, conversely, at a critical interaction 
U C 2{T)). This third phase turns out to be less stable as 
compared to the metal and to the insulator in the entire 
parameter regime. Similarly, one can define a critical 
temperature T C \{U) (at smaller U) and thus a critical 
interaction U C \{T) where the insulating phase coalesces 
with the third phase. 

Calculations for different U and T have been performed 
to obtain the entire phase diagram within the two-site ap- 
proach. The result is shown in Fig. [3] For U < U C 2(T) 
there is a metallic phase which is smoothly connected to 
the U = limit. On the other hand, for U > U cl {T) there 
is an insulating phase which is smoothly connected with 
the Mott insulator for U — ► oo. Metallic and insulating 
phase are coexisting for U C \(T) < U < U C 2(T). At zero 
temperature, the metal is stable as compared to the insu- 



lator in the entire coexistence region, and the transition 
is continuous at U c = U C 2- For finite temperatures the 
transition is discontinuous at a critical interaction U C (T) 
with U cl (T) < U C (T) < U c2 (T). With increasing tem- 
peratures the coexistence region U C \(T) < U < U C 2(T) 
shrinks and disappears above a critical temperature T c 
which is defined via U cl (T c ) = U c2 (T c ). Above T c the 
metallic phase is smoothly connected with the insulating 
phase. 

Qualitatively, this is exactly the same result that is ob- 
tained within the full DMFT. 11 It is very remarkable that 
the rather complex topology of the phase diagram can be 
reproduced with a comparatively simple two-site model 
as a reference system. This shows that it is essential to 
perform the mapping onto the reference system in a way 
which is thermodynamically consistent and which is con- 
trolled by a physical variational principle. Note that the 
L-DMFT, or its extension away from the critical point 
at T = 0^ fails to reproduce the discontinuous transi- 
tion for T > and the critical temperature T c due to the 
ad-hoc character of the approximation. 

Quantitatively, one should expect some deviations 
from the results of the full DMFT due to the simplicity 
of the two-site model. Comparing with the NRG result^ 
for U C \(T = 0) « 4.8 — 5.0, the two-site approximation 
overestimates the critical interaction by a few per cent. 
The determination of the critical temperature is difficult 
in any numerical approach. T c at 0.05 — 0.08 is estimated 
from the QMC and NRG results of Refs. l2l27i Thus, the 
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two-site approximation underestimates T c by more than 
a factor 2. It is worth mentioning that the numerical ef- 
fort to obtain the entire phase diagram in the U—T plane 
is of the order of a few minutes on a standard worksta- 
tion which is negligible as compared to a DMFT-QMC 
calculation. 



IX. CONCLUSIONS 

Our present understanding of the Mott-Hubbard tran- 
sition is mainly based on the exact solution of the one- 
dimensional case on the one hand and on the D = oo 
mean-field picture provided by the DMFT on the other 
hand. For the physically more relevant two- or three- 
dimensional Hubbard model, however, neither the ana- 
lytical concepts developed for D = 1 nor the mean-field 
theory can be expected to give an essentially correct and 
comprehensive description. Direct numerical approaches, 
such as QMC, are able to give essentially exact results 
for a D — 2, 3 dimensional lattice of finite spatial extent. 
However, the relevance of the results for the thermody- 
namic limit and, in many cases, for the low-temperature, 
low-energy regime remains an open question. 

In this situation, a combination of a direct numerical 
approach for systems of finite size with the mean-field 
concept appears to be advantageous. This is more or less 
the direction that is followed up by the different clus- 
ter extensions of the DMFTi 37 i 51 i 52 i 53 Via a generalized 
mean- field concept, the original lattice problem given by 
the Hamiltonian H is mapped onto a cluster problem 
described by a Hamiltonian H' . In fact, due to the pres- 
ence of strong short-range antiferromagnetic correlations 
a considerable revision of the mean-field picture is prob- 
ably necessary^ However, the cluster extensions of the 
DMFT suffer from the fact that the mean-field formula- 
tion requires that the sites of the finite cluster are coupled 
to uncorrelated baths with an infinite number of degrees 
of freedom each. This circumstance complicates the prac- 
tical solution of the H' problem (which must be solved re- 
peatedly during the self-consistency cycle) so much that 
additional approximations are required and/or stochastic 
numerical methods. 

The self-energy-functional approach offers an interest- 
ing alternative as the reference system H', the original 
model H is mapped onto, is by no means completely pre- 
determined. The SFA gives a very general prescription 
how this mapping can be performed while keeping the 
thermo dynamical consistency of the approach as an ex- 
plicit expression for the grand potential is provided. In 
this way the DMFT and the cellular DMFT are recov- 
ered as certain limits for special choices of H' , namely for 
a decoupled system of clusters of size iV c = 1 or N c > 1 
including couplings to rib = oo bath degrees of freedom. 
There is, however, the additional possibility to construct 
approximations with rib < oo which are consistent in 
themselves in the same way as are the DMFT and the 
C-DMFT. 



Now, the question is whether or not an approximation 
with a finite number of bath sites rib < oo can be toler- 
ated. Note that even in the (C-)DMFT it becomes nec- 
essary to reduce the problem posed by W with rib = oo 
to a numerically tractable one with a finite number of 
degrees of freedom. While this additional approxima- 
tion is controlled within QMC or ED approaches, for ex- 
ample, it nevertheless violates thermodynamical consis- 
tency. Within the SFA, on the contrary, the approxima- 
tion is derived from a thermodynamical potential for any 
rib, infinite or finite or even so small as rib = 1. Depend- 
ing on the quantity and accuracy one is interested in, 
there can be a rapid convergence with respect to rib (cf. 
Ref. l3Cf) so that a small number might be sufficient. For 
cluster approximations, the best choice is by no means 
clear as it must be balanced with the choice of the cluster 
size N c . This strongly depends on the lattice dimension. 
In fact, it has been shown 3 ^ that for the one-dimensional 
Hubbard model a larger N c is to be preferred as com- 
pared to a larger rib- 

In this context, it is interesting to see where one are 
led to with the most simple reference system conceiv- 
able. This is a model H' characterized by N c = 1 and 
?ib = 1 which yields the so-called two-site dynamical im- 
purity approximation (DIA). The answer is given with 
the present paper: Even in this approximation the Mott 
transition shows up. At zero temperature the transition 
turns out to be continuous at a finite critical interac- 
tion U c where there is a coalescence of the metallic with 
a coexisting insulating phase. For finite temperatures, 
on the other hand, the transition is discontinuous. The 
first-order line U C (T) terminates at a second-order criti- 
cal point (U C (T C ),T C ), and for T > T c there is a smooth 
crossover only. This is qualitatively the same picture 
that has been found beforehand in the full DMFT. Fur- 
thermore, the two-site DIA yields a U c at T = that 
is surprisingly close to the result of the full DMFT. In 
this respect the two-site DIA is of similar quality as the 
linearized DMFT^ another approach that is based on 
a mapping onto the two-site SIAM. Whether or not the 
t/ c -results of the two-site DIA improve on those of the L- 
DMFT is difficult to decide in view of the existing (full) 
DMFT data. More important, however, the conceptual 
improvement gained is substantial: While the mapping 
procedure on the two-site model, though physically mo- 
tivated, is done in an ad-hoc way in the L-DMFT, the 
two-site DIA is derivable from a thermodynamical po- 
tential and can be characterized as an optimal two-site 
approach in fact. 

The fact that a reasonable mean-field description of 
the Mott transition is possible even for the simple two- 
site DIA, motivates further SFA studies of the transi- 
tion with improved approximations in the future (e.g. 
using finite clusters, N c > 1, and small rib). In case 
of larger and more complex reference systems, however, 
it becomes more and more important to have a numeri- 
cally efficient and accurate method for the evaluation of 
the self-energy functional at hand. The detailed analysis 
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of the different contributions to the functional has been 
another major intention of the present study. Using the 
causality properties of the Green function, it has been 
shown that the important Trln(— G) term can be writ- 
ten as the grand potential of a system of non-interacting 
quasi-particles with unit weight (apart from a correc- 
tion term that cancels out in the functional eventually). 
While this result is also interesting by itself, it is very 
well suited for an efficient numerical implementation. As 
the energies of the ficticious quasi-particles are given by 
the poles of the Green function, they can be found in a 
comparatively simple way by exploiting general causality 
properties. In particular, there is no need for a small 
but finite Lorentzian broadening, lo — > lo + iS, which 
must be introduced in a more direct way to evaluate the 
functional This will become important when study- 
ing the critical regime of the Mott transition, where an 
accurate computation of small energy differences is vital, 
using an approximation with several variational parame- 
ters. Studies in this direction arc intended for the future. 
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with C, D Hcrmitian and D = B 2 positive definite. 
Let C = UcU^ with U unitary and c diagonal. Then 

-. ] — = DU —J- U f D = DU 4^r = X ^iY 

A±iB c±i c 2 + 1 

with X Hcrmitian and 



Y = DU 



U f D 



(A3) 
(A4) 



Hermitian and positive definite. 

Consider a fixed (real) frequency lo. Since Go is causal, 
G ,rct = G ,R-iG ,i with G ,r, G ,i Hermitian and G ,i 
positive definite. Using the lemma, Gq ret = Pr + iP\ 
with Pr , Pi Hcrmitian and P\ positive definite. Since S 
is causal, S rc t = Sr — i^i with Sr, Si Hermitian and 
Si positive definite. Therefore, 



G rc 



1 



1 



Pr, + iPi - S R + iSi Q R + iQi 



(A5) 



with Qr Hcrmitian and Qi Hermitian and positive defi- 
nite. Using the lemma once more, shows G to be causal. 
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APPENDIX A 

Here is will be shown that the Green's function G is 
causal: For any t' the self-energy X = E(t') is causal 
since it is defined to be the exact self-energy of the 
reference system H' = Ho{t') + H\(U). Likewise, the 
Green's function Go is causal. It has to be shown that 
G = (Gq -1 — S) _1 is causal if £ and the free Green's 
function Go are known to be causal. 

Causality of G means (i) that G a /3(u)) is analytic in 
the entire complex lo plane except for the real axis and 
(ii) that G re t(w) = G(lo + i0 + ) = Gr, — iG\ for real lo 
with Gr, Gi Hcrmitian and Gi positive definite, (i) is 
easily verified. To show (ii) we need the following 

Lemma: For Hermitian matrices A, B with B positive 
definite, one has 



1 



A±iB 



(Al) 



with X , Y Hcrmitian and Y positive definite. The proof 
of the Lemma is straightforward: 



1 



A±iB 



= B "4 



1 



BiABi±i 



B 5 = D 



1 



C±i 



D 

(A2) 



APPENDIX B 



Consider the function 



with real isolated first-order poles at lo 
residues R m > 0. We have: 



(Bl) 



and positive 



-Im/(w + i0 H 

7T 



> 



(B2) 



In Sec-EH/M = g k (uo) = 1/(lo + »- Vk (u)). 

In a polar representation —f(u> + i0 + ) = r{Lo)e 1 ^^ 
with — 7r < 4>{lo) < 7r. On the principal branch of the 
logarithm and for real lo one has Imln(— /(lo + i0 + )) = 
4>(lo). For any lo ^ Lo m we have Im(— f(u> + «0 + )) = + 
from Eqs. (jBip and l|B2|l and thus <)){uj) = ir for —f(oo) < 
0, and 4>(lo) = for —/(lo) > 0. Consequently, 

Imln(-/( W + i0+)) = 7T 0(/(w)) = 7T 0(l//(w)) , (B3) 

where Q is the step function. For lo = cj m , the imaginary 
part of — f(u> + i0 + ) diverges. However, — ir < 4>(lo) < ir 
(in fact (j)(Lo m + 0+) = -vr/2 and ^(w m - 0+) = tt/2). 
Hence, Im ln(— f(u) + i0 + )) remains finite at lo = Lo m and 
can be ignored in an integration over real lo as the poles 
of f(co) are isolated. 



APPENDIX C 

Consider a function / = /(lo) which is analytical ex- 
cept for isolated first-order poles on the real axis and 
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which is real for real u>. 0(/(w)) is constant almost ev- 
erywhere, and thus the derivative d@(f(u>))/du> vanishes 
almost everywhere. A non-zero derivative may either oc- 
cur if f(u>) = 0, i. e. at the zeros u) m . This gives a 
contribution: 



<5(/M)J» 



? I/'KJI 



S(u> - U} m ) 



(CI) 



where f'(u>) = df(uS)/duj. A second contribution arises 
from the first-order poles of f(ui) at £„. The poles are 
the zeros of the function l/f(u>). Note that 0(/(cj)) = 
0(1// (w)) since the sign of the argument is unchanged. 
Thus the contribution due to the poles is: 

<j(i//M)(i//)'M = £ - Cn) (C2) 



We have: 

<*e(/M) 
do; 



<5(/H)/'H + < 5(i//M)(i//)'^) 

{zeros} {poles} 

e m 5(u; - cj m ) + ^ e n 8(uj - Cn) 

m n 

(C3) 



with e m = ±1 and e n — ±1 depending on the sign of the 
slope of / at the zeros uj m and the sign of the residue of 
/ at the poles ( n , respectively. 

Consider the function f(to) = lu + [i — rik{<jS). The 
zeros of / are the poles of the diagonalizcd Green's func- 
tion l/(w + /U — r?fc(w)) which has positive residues. This 
implies a positive /' at oj m and e m = +1. The poles 
of / are the poles of the self-energy which has positive 
residues. Thus the residues of / at the poles are negative 
and (1//)' is negative at and e„ = — 1. Hence: 



du 



£%-c«). 

n 

(C4) 



APPENDIX D 

The reference system H' is a set of decoupled single- 
impurity Anderson models with one-particle energy of 
the impurity site e\, conduction-band energies €k with 
k = 2, ...,n s and hybridization strengths Vk- We first 
calculate the eigenvalues of the hopping matrix t' . The 
matrix is block diagonal with respect to the site index i. 
Each block is labeled by the "orbital" index k = 1, ...,n s . 
There are non-zero elements of the matrix for k = 1, 
k' = 1, and k = k' only: 



Furthermore, Yi^{uj) 




5fci<5/c'iS 




Using 








/ di a\ a\ a\ . 


;^ 




/ di 












a 2 d 2 . 









d 2 







det 


a 3 d 3 . 




= det 





d 3 









04 di . 












di 






v 






u 









(D2) 

where di = di — X)fe=2 \ a k\ 2 /dk and the general relation 
trlnA = In det A, we find: 

itr'ln(-G'M) = In det ~* — r 



In det 



/ -Gi 

-G' 2 







-G' 



\ 



ln(-G 1 M)+^ln(-G fe M) 

(D3) 



with G\{iuj) and Gk(ioj) for fc = 2,...,n s as defined by 
Eqs. lfTT|l and l(T2*)l. The factor 1/2 accounts for the two 
spin directions. 
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